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Abstract: Osteosarcoma (OS) is the most common primary malignant tumor of the bone, 
and pulmonary metastasis is the most frequent cause of OS mortality. The aim of this study 
was to discover and characterize genetic networks differentially expressed in metastatic 
OS. Expression profiling of OS tumors, and subsequent supervised network analysis, was 
performed to discover genetic networks differentially activated or organized in metastatic 
OS compared to localized OS. Broad trends among the profiles of metastatic tumors include 
aberrant activity of intracellular organization and translation networks, as well as 
disorganization of metabolic networks. The differentially activated PRKC8-RASGRP3-GNB2 
network, which interacts with the disorganized DLG2 hub, was also found to be 
differentially expressed among OS cell lines with differing metastatic capacity in xenograft 
models. PRKCe transcript was more abundant in some metastatic OS tumors; however the 
difference was not significant overall. In functional studies, PRKCe was not found to be 
involved in migration of Ml 32 OS cells, but its protein expression was induced in Ml 12 
OS cells following IGF-1 stimulation. 
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1. Introduction 

Osteosarcoma (OS) is the most common primary malignant cancer of the bone, with nearly 
1,000 diagnoses annually in North America [1,2]. High-grade intramedullary osteosarcoma 
("conventional") constitutes the majority of cases, and is an aggressive disease that typically 
metastasizes to the lungs [1,2]. Approximately 10-20% of patients present to the clinic with 
discernable metastasis, and a further 20-30% develop metastasis despite aggressive treatment [3]. 
Currently neoadjuvant chemotherapy in combination with surgical resection of osteosarcomas can 
achieve 5-year overall survival rates of approximately 65% [4,5]. However, delineation of the molecular 
mechanisms contributing to osteosarcoma metastasis has the possibility of identifying specific 
therapeutic targets that improve prognosis of patients with metastatic osteosarcoma. 

Expression profiling is an "omic"-level discovery technique that allows the quantification of 
thousands of different transcripts simultaneously between different disease states [6]. This technology 
has many applications in the cancer field including distinguishing subtypes of a particular cancer, 
identifying transcripts differentially expressed between those types, and predicting the subtype of a 
cancer based on the expression of identified transcripts [6,7]. In recent years it has been shown that 
analysis of expression profiles at the level of sets of transcripts rather than individual transcripts has 
more statistical power, is more reproducible between studies, and provides informative context to the 
results [7-9]. There are many algorithms for gene-set or network-based analytical algorithms, 
including the popular gene-set enrichment analysis (GSEA) and the proprietary ingenuity pathway 
analysis (IP A) [10]. Additionally dozens of other network-based analytical algorithms have been 
published that either utilize different statistical methods or that interrogate a unique cellular application 
(reviews: [8,9]). Two such algorithms are the methods of Chuang et al that aims to discover 
differentially "activated" genetic networks [11], and "Dynemo" published by Taylor et al. that aims to 
discover differentially organized genetic networks [12]. These algorithms have successfully been 
applied to uncover molecular alterations in expression profiles of poor-outcome breast cancers, and 
have achieved over 70% predictive accuracy in -classifying those tumor samples utilizing a five-fold 
cross-validation strategy [11,12]. 

The aim of the present study was to identify differentially activated and organized networks in 
expression profiles of metastatic-at-diagnosis osteosarcomas (MD-OS) compared to expression profiles 
of localized-at-diagnosis osteosarcomas (LD-OS). This analysis would assist in the prioritization of 
candidate networks for in vitro characterization in osteosarcoma cell lines, including cells of differing 
metastatic capability. 
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2. Results 



2.7. Unsupervised Hierarchical Clustering of Expression Profiles Reveals Distinct Subtypes of 
Osteosarcomas 

Expression profiling was performed on 46 LD-OS tumor samples and 17 MD-OS tumor samples 
using UHN "human 19 k cDNA microarrays (H19K)", that contain approximately 19,000 cDNA spots 
mapped to approximately 9,000 unique transcripts. Unsupervised gene discovery techniques were used 
to identify a subset of 596 genes that exhibited at-least six-fold change in expression in at-least four 
tumors (p < 0.001), and hierarchical clustering of tumors according to the expression pattern of these 
genes demonstrated two broad groups of tumors. One group (Figure IB, red dendrogram) contained all 
but two of the MD-OS tumors, and the other group (Figure IB, green and blue dendrograms) contained 
all but one of the LD-OS tumors. This analysis supports the notion that osteosarcomas of varying 
metastatic status at the time of diagnosis exhibit distinct expression patterns, and that division of 
tumors into these categories (localized at diagnosis or metastatic at diagnosis) is appropriate for 
subsequent supervised analysis. 

Figure I. (A) Unsupervised hierarchical clustering of 596 highly variable genes. 
(B) Dendrogram depicting separation of MD-OS tumors (green bars) and LD-OS tumors 
based on expression of highly variable genes. (A) Unsupervised hierarchical clustering of 

05 tumor samples (Partek Genomics Suite) was performed on a subset of 596 genes found 
to be significantly differentially expressed within the dataset (expression varies by at least 

6 fold, in at least 4 patients, p < 0.001), according to methods previously described [13,14]. 
(B) The dendrogram reveals two clusters, one with red dendrogram lines that contains 
mostly metastatic-at-diagnosis tumors (green rectangles). A second cluster, which can be 
subdivided further (blue and green dendrogram lines) contains mostly localized-at-diagnosis 
tumors (purple rectangles). 



A 




B 



lllllllllllllllllllllllllllllllllllllllllllllllllllllllllllllll 



Cancers 2013, 5 



375 



2.2. Differentially Activated and Organized Networks in Metastatic Osteosarcomas 

The expression information was integrated with the Pathway Commons database of physical and 
genetic interactions [15], yielding a dataset of 5,855 genetic networks. Supervised analysis was 
performed on these data in a manner similar to previous publications to identify the differentially 
activated and organized networks in MD-OS expression profiles. To describe broad trends in human 
metastatic OS, inclusive cut-off levels were identified to delineate as many networks as possible that 
were significantly differentially activated or organized in the metastatic samples. For this analysis, four 
hundred and ninety seven (497) of the 5,855 genetic networks were identified as differentially 
activated and six hundred and eighty three (683) networks were significantly differentially organized. 
Networks annotated to transport, translation, organization and protein modification were more 
commonly differentially activated than expected by chance (Table 1, column 8). Networks annotated 
to the processes of organization, transport and translation were more commonly differentially organized 
than expected by chance, as were metabolic networks (Table 1, column 11). Three hundred thirty eight 
(338) differentially activated, and one hundred and sixty two (162) differentially organized, networks 
were visualized and clustered according to their Gene Ontology Process annotation (Figure 2). This 
visualization allowed the observation that clusters of differentially activated networks interact with large, 
significantly disorganized hubs (e.g., Figure 2 — "protein modification" process: cluster of 
differentially activated networks interacting with the disorganized hub breakpoint cluster region 
(BCR). This pattern is also observed in the "translation" and "transport" processes. 

Table 1. Enrichment of genetic networks in metastatic osteosarcomas according to cellular 
processes at the global (Permissive) significance level. 

Cellular Processes All Networks in Study Significant Networks: Permissive Cutoffs 



Subsets Entire Networks Differentially Activated Differentially Organized 



1 


2 


3 


4 


5 


6 


7 


8 


9 


10 


11 




# 


% of all # 


% of all 


# 


% of all 


p -value 


# 


% of all 


p -value 


protein modification 


573 


9.8 


575 


9.8 


52 


10.5 * 


5.3E-02 


49 


7.2 




transport 


1002 


17.1 


869 


14.8 


194 


39* 


l.lE-34 r 


133 


19.5 * 


9.0E-05 r 


signaling 


954 


16.3 


955 


16.3 


69 


13.9 




87 


12.7 




transcription 


595 


10.2 


635 


10.8 


17 


3.4 




48 


7 




stress 


247 


4.2 


257 


4.4 


15 


3 




9 


1.3 




metabolism 


1340 


22.9 


1410 


24.1 


58 


11.7 




275 


40.3 * 


5.7E-24 r 


cell cycle 


65 


1.1 


68 


1.2 


6 


1.2 * 


1.6E-01 


7 


1 




reproduction 


130 


2.2 


137 


2.3 


6 


1.2 




0 


0 




intracellular organization 


195 


3.3 


201 


3.4 


27 


5.4* 


3.4E-03 r 


54 


7.9 * 


9.1E-10 r 


development 


338 


5.8 


357 


6.1 


0 


0 




0 


0 




translation 


146 


2.5 


99 


1.7 


50 


10.1 * 


4.4E-19 r 


17 


2.5 * 


3.9E-2 r 


death 


31 


0.5 


36 


0.6 


0 


0 




3 


0.4 




cytoskeleton 


57 


1 


59 


1 


0 


0 




0 


0 




ion transport 


10 


0.2 


8 


0.1 


2 


0.4* 


1.6E-01 


0 


0 




proliferation 


3 


0.1 


2 


0 


1 


0.2 * 


2.1E-01 


0 


0 




homeostasis 


2 


0 


2 


0 


0 


0 




0 


0 




differentiation 


6 


0.1 


6 


0.1 


0 


0 




0 


0 




NaN 


161 


2.7 


179 


3.1 


0 


0 




1 


0.1 




Totals 


5855 


100 


5855 


99.8 


497 


100 




683 


99.9 





* denotes an enrichment (increase in proportion) above that found in the appropriate "all networks in study" 
category (e.g., column 7 is compared to column 3, and column 10 is compared to column 5). r denotes a 
significant /?-value (p < 0.05). 
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Figure 2. Nodes represent differentially activated (red squares) and organized (blue 
squares) networks in metastatic OS samples. Green edges represent genetic overlap 
between networks. Genes in significant networks were annotated with simplified Gene 
Ontology Slim Generic terms (Table 6), and networks were grouped into processes by the 
most commonly occurring term in each network. Networks meeting the cut-off conditions 
detailed at the top of the figure were visualized with the Enrichment Map plugin for 
Cytoscape. The size of each node reflects the number of genes included in the network. 




Node size: network size Edge width: degree of overlap 6 genes 100 genes 500 genes 

2.3. Genes Previously Implicated in Osteosarcoma Metastasis are among Significant Network Results 

Next, it was examined if the genetic networks discovered to be differentially activated or organized 
in expression profiles of MD-OS samples contained genes previously implicated in OS metastasis. 
Extensive literature review led to the generation of a query list containing genes whose expression is 
either correlated with outcome in human OS patients, or genes that have been shown to be involved in 
OS metastasis in xenograft studies or animal models of OS (Table 2). Thirty-two genes from this list 
were included in the present study (i.e., present on the expression profiling microarrays), and it was 
found that five such genes were among networks found to be significantly differentially activated in 
MD-OS expression profiles, and twenty-eight were among networks found to be significantly differentially 
organized in MD-OS expression profiles (Tables 3 and 4). 



Table 2. Genes previously implicated in metastatic progression of OS. 



Gene 


Ref. 


Gene 


Ref. 


Gene 




Ref. 


BIRC5 


[16] 


IVD 




CXCR3 


X 


[17] 


CAV1 


[18] 


KIT 


[19] 


EGF 


X 


[20,21] 


CCN1 


[22] 


LRP5 


[23] 


EGFR2 


X 


[24-29] 


CCN3 


[30] 


MMP2 


[31] 


FADD 


X 


[32] 


CD44 


[33] 


PDGFRA 


[34,35] 


HIF1 


X 


[36] 
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Table 2. Cont. 



Gene 


Ket. 


Gene 


Ket. 


Gene 




Ket. 


CDH11 


[37] 


PEDF 


ro on 

[38] 


TT 1 ^> A 

1L12A 


X 


[39] 


CDH2 


[37] 


RECK 


[40] 


IL12B 


X 


[39] 


COL18A1 


r a 1 n 

[41] 


S100A6 


r a 

[42] 


INS 


X 


r a o c 1 n 

[43-51] 


CXCR4 


re ^ r n 

[52-55] 


SPARC 


[56] 


A /f A A 1 

MAML1 


A 


[57] 


DPF2 




TIMPl 


[58] 


MIRK 


X 


[59] 


EGFR 


[60] 


PLAUR 


[61] 


MMP14 


X 


[62] 


EZR 


[25,63,64] 


VEGFA 


[55,65-72] 


MMP9 


X 


[58,73,74] 


FAS 


[32,75-77] 


VEGFB 


[55,65-72] 


PDGFA 


X 


[35] 


IGF1 


[43-51] 


VEGFC 


[55,65-72] 


VCP 


X 


[78] 


IGF2 


[43-51] 


WIFl 


[79] 








IGF1R 


[43-51] 


YYI 


[80] 









X denotes genes not present in the expression profiling microarrays used in this study. 



Table 3. Differentially activated networks containing genes previously implicated in 
OS metastasis. 



Query 


Node 


Score 


Psample 


FDRsample 


Pgene 


FDRgene 


CDH2 


CDH2 


3.89 


0 


0.0017 


0.037 


0.019 


S100A6 


CACYBP 


3.94 


0 


0.00096 


0.007 


0.18 


TIMPl 


ECH1 


3.28 


0 


0.00096 


0.003 


0.17 


PLAUR 


PGAP1 


2.54 


0 


0.0019 


0.003 


0.19 


IVD 


MECR 


3.94 


0.001 


0.01 


0.046 


0.0076 



"Node" is the central gene exhibiting differential network activity. 



Table 4. Differentially organized networks containing genes previously implicated in 
OS metastasis. 



Query 


Node 


APCCxotal 


Psample 


Query 


Node 


APCCxotal 


Psample 


S100A6 


S100A6 


0.41 


0 


TIMPl 


LRP1 


0.32 


0 


KIT 


KIT 


0.34 


0 


PLAUR 


LRP1 


0.32 


0 


YYI 


YYI 


0.32 


0.001 


BIRC5 


PAFAH1B1 


0.34 


0 


BIRC5 


BIRC5 


0.35 


0.001 


NOV 


GIA1 


0.3 


0 


MMP2 


MMP2 


0.35 


0.001 


VEGFA 


SPARC 


0.34 


0 


SPARC 


SPARC 


0.34 


0 


YYI 


HDAC 


0.33 


0 


IGF1R 


CAMK2B 


0.37 


0 


IGF1 


PRKCD 


0.36 


0 


IGF2 


CAMK2B 


0.37 


0 


PDGFRA 


JAK1 


0.33 


0 


EGFR 


CAMK2B 


0.37 


0 


CXCR4 


JAK1 


0.33 


0 


KIT 


CAMK2B 


0.37 


0 


FAS 


BTK 


0.35 


0 


CDH2 


CAMK2B 


0.37 


0 


VEGFB 


RASA1 


0.32 


0 


EZR 


CAMK2B 


0.37 


0 


MMP2 


ITGB2 


0.37 


0 


CD44 


CAMK2B 


0.37 


0 


LRP5 


FZD8 


0.48 


0 


SPARC 


CAMK2B 


0.37 


0 


SERPINF1 


CSNK2A1 


0.33 


0.001 


CYR61 


ATP2A2 


0.33 


0 


COL18A1 


CTSL1 


0.33 


0.001 


CAV1 


ATP2A2 


0.33 


0 


WIFl 


FZD1 


0.37 


0.001 


S100A6 


ACTN1 


0.32 


0 


IVD 


MCCC1 


0.37 


0 



"Node" is the central gene exhibiting differential network organization. 
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2.4. The PRKCs-RASGRP3-GNB2 Network Is Differentially Activated, and May Interact with the 
Disorganized DLG2 Hub 

In order to identify drivers of the metastatic phenotype, the lists of significant networks discovered 
to be differentially activated and organized in expression profiles of MD-OS samples were further 
refined by focusing on a manageable number of each type of network that exhibited the greatest 
change between LD-OS and MD-OS samples. Twelve differentially organized networks were identified 
that exhibited the greatest change in organization between LD-OS and MD-OS samples, and that were 
also "genetic hubs" (having more interactors than the median), as it has been shown that these genes 
may be particularly important to cellular processes (Figure 3) [81]. Forty-three networks were 
identified that exhibited the greatest change in activity between LD-OS and MD-OS expression 
profiles, and their interactions with some of the most disorganized hubs were visualized (Figure 4). 
This analysis allowed the identification of the networks exhibiting the greatest change among MD-OS 
expression profiles, and also the observation that many differentially activated networks also interact 
with significantly disorganized hubs, e.g., the differentially activated network containing protein 
kinase C epsilon (PRKCe), RAS guanyl releasing protein 3 (RASGRP3), and guanine nucleotide 
binding protein 2 (GNB2), also interacts with the disorganized network of discs large homolog 2 
(DLG2) (Figure 4). 

Figure 3. Labeled black squares show 12 differentially organized hubs (APCC > 0.4, p = 0, 
interactors > 7). Nodes are genes and their colour reflects fold change in expression in the 
metastatic samples, only for significantly differentially expressed genes (t-test p < 0.001). 
Only interactors of a significant hub whose expression is significantly correlated (PCC 
p < 0.001) with expression of that hub in either the localized or metastatic samples are 
included. Edges depict interactions, and their colour reflects change in co-expression 
between metastatic samples and localized samples (green: correlation increased in 
metastatic samples, purple: correlation decreased in metastatic samples). Asterisks denote 
co-expressions that are significant (PCC p < 0.001) in both localized and metastatic samples. 




Nodes: genes Edges: : interactions that are 

^^^^^^^^^^ ^^^^^^^^^^ significant in both 

LooMets Loc < Mets LooMets Loc < Mets 
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Figure 4. Genes comprising 43 differentially activated networks that meet stringent cut-offs 
(detailed at the top of the figure) were visualized with Cytoscape, as well as their 
interactions with 11 differentially organized networks (that also met stringent cut-offs 
which are detailed at the top of the figure). Coloured nodes represent relative abundance in 
metastatic samples vs. localized samples (red = more abundant, blue = less abundant). 
Coloured edges represent correlations between genes in metastatic samples vs. localized 
samples (green = more correlated, purple = less correlated). The number of included interactors 
of each gene in this study is proportional to the size of the node depicting the gene. 



•43 most differentially active Chuang Networks 
[ J 11 most differentially organized Dynemo Hubs 



P's<0.01, FDR's< 0.1, ANAS >1.65 
P = 0, APCC ^0.4, interactors:: 6 U sources) 



IKBKE 



MCC 



HIPl 



DLG2 



PRKCE 

RASGRP3 

GNB2 




Nodes: genes 
Colour: expression 



n 

Edges: interactions 
Colour: correlation 



ooMets 



Loc < Mets Loc > Mets 



Loc < Mets 



Representative Sizes 
10 interactors 500 interactors 



Literature searches of the genes contained in the highly significant network results showed that 
PRKCe has been characterized in progression of prostate [82-87], breast [88,89], and renal cancers [90], as 
well as tumorigenesis of squamous cell carcinoma (non-melanoma skin [91-95], and head and 
neck [96-99]), as well as non-small cell lung cancer [100-103]. Additionally high expression of GNB2 
is associated with an aggressive form of pulmonary adenocarcinoma (mixed adenocarcinoma with 
bronchioalveolar features), and shorter overall survival for patients with these tumors [104]. RASGRP3 
has been shown to promote androgen independence and progression of prostate cancer [105]. These 
reports of the involvement of these genes in the progression of other cancers led to the selection of the 
differentially activated PRKCs-RASGRP3-GNB2 network as a lead candidate for further characterization 
in osteosarcoma cells and tumors (Figure 5). 
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Figure 5. (A) The PRKCe-RAS GRP3 -GNB 2 Network is Differentially Activated in 
Expression Profiles of Metastatic Osteosarcomas. (B) Differential Activity Score of 
PRKC8-RASGRP3-GNB2 Network. (A) Shown is the network of PRKCs-RASGRP3-GNB2 
that is differentially activated in metastatic OS samples of the expression profiling screen. 
Red nodes are transcripts that are more abundant in metastatic samples, and blue 
nodes are less abundant in metastatic samples. Within each circle is the %-fold change in 
expression of metastatic samples vs. localized samples, and the order of magnitude of the 
p-value from the corresponding student's t-test. (B) Shown are the relevant scores for the 
PRKC8-RASGRP3-GNB2 network, the results of the differential network activity analysis 
applied to expression profiles of MD-OS vs. LD-OS tumors. 




Network PI FORI P2 FDR2 

Activity 

(sample (sample (gene (gene 

core (NAS) 

permutations) permutations) permutations) permutations) 



Value 3.31 0 0.0011 0.008 0.0878 



2.5. PRKCs-RASGRP3-GNB2 Network Is Differentially Activated in Vitro 

A panel of human osteosarcoma cell lines known to have differing metastatic potential when grown 
as murine xenografts was collected to investigate the role of PRKCs-RASGRP3-GNB2. This panel 
includes the parental Hu-09 cell line (included in Figure 6), and its highly metastatic derivatives Ml 12 
and M132 (derived by in vivo metastatic selection), as well as the poorly metastatic sub-clones L06 and 
L13 (all generously provided by Dr. B. Fuchs, the University of Zurich, Zurich, Switzerland) [106,107]. 
Additionally the poorly metastatic SAOS2 and MG63 cell lines, as well as their highly metastatic 
derivatives LM7 and M8, were included (generously provided by Dr. E. Kleinerman, University of 
Texas MD Anderson Cancer Center, Houston, TX, USA) [108,109]. In this panel, it was observed that 
the PRKC8-RASGRP3-GNB2 network exhibited an mRNA expression pattern similar to that observed 
in the expression profiles of osteosarcoma tumors (Figure 6). Specifically, it was observed that PRKCs 
mRNA was significantly more abundant in some highly metastatic lines (Ml 12 and Ml 32 vs. L06 and 
LI 3), and that RASGRP3 was significantly less abundant in some highly metastatic lines (LM7 vs. 
SAOS2 and M8 vs. MG63). Additionally, mRNA levels of myosin chain heavy 9 (MYH9), which has 
been shown to interact with PRKCe at stress fibers in mice [110], was also observed to be more 
abundant in some highly metastatic lines (LM7 vs. SAOS2 and M8 vs. MG63). The protein expression 
of PRKCe was also investigated in this panel, and it was found that the highly metastatic lines Ml 12 
and Ml 32 expressed PRKCs protein at higher levels than observed in either the less metastatic 
parental Hu09 line, or in the poorly metastatic LI 3 derivative (Figure 7). 
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Figure 6. The PRKC8-RASGRP3-GNB2 network is differentially activated in vitro. 
Quantitative PCR was performed on cDNA synthesized from the human OS cell lines 
shown. Purple bars denote cell lines that are poorly metastatic in mouse models, and green 
bars denote cell lines that are strongly metastatic as murine xenografts. Arrows denote 
relationships between cell lines. Student's t-test: * p < 0.05 ** p < 0.001. Multiple 
independent experiments are shown for each gene: GNB2 (n = 2), PRKCs and RASGRP3 
(n = 3), and MYH9 (n = 5). ** p < 0.001, * p < 0.05. 
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Eft: 
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Figure 7. PRKCs Protein is more abundant in some highly metastatic cell line models of 
osteosarcoma metastasis. Western blots of whole cell ly sates from human OS cell lines 
probed with anti-PRKCs and anti-(3-actin antibodies. Purple bars denote cell lines that are 
poorly metastatic in mouse models, and green bars denote cell lines that are strongly 
metastatic in cell lines. Arrows denote relationships between cell lines. Shown is a 
representative western blot of two independent experiments. 
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2.6. Human Osteosarcomas That Are Metastatic-at-Diagnosis Are More Likely to Exhibit High Levels 
of PRKCs mRNA 

The amount of PRKCs mRNA was assayed in 17 LD-OS and 14 MD-OS tumors that were used in 
the expression profiling screen (Figure 8). It was found that MD-OS tumors do not have a higher 
average expression of PRKCs than LD-OS tumors (Welch two sample t-test p = 0.1873), however 
MD-OS tumors were more likely to exhibit high expression of the PRKCs transcript (five of fourteen 
MD-OS tumors: Figure 8). 

Figure 8. Metastatic-at-Diagnosis human osteosarcomas are more likely to have high 
PRKCs expression. Quantitative PCR was performed on cDNA synthesized from tumor 
samples from the original expression profiling cohort. PRKCs expression was normalized to 
STAM2, and is depicted as fold-change. The purple bar denotes localized tumors at 
diagnosis (LD-OS), and the green bar denotes metastatic tumors at diagnosis (MD-OS). 
The dotted red line depicts the average PRKCs expression in localized tumors. The number 
of replicates for each sample is shown below each bar. 
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2. 7. PRKCs Is Not Required for Migration of Highly Metastatic Ml 32 Cells 

Since PRKCe promotes the in vitro migration of other cell lines [89], and the in vitro migration rate 
of the Hu09-derived cell lines correlates with their ability to form metastatic colonies in mice [107], 
experiments were undertaken to determine whether PRKCe promotes migration of osteosarcoma cells. 
The highly metastatic Ml 32 cell line was selected for PRKCe-knockdown studies, as this line 
expressed PRKCs protein at a high level (Figure 7). However, it was observed that knockdown of 
PRKCs protein using siRNA did not affect the in vitro migration rates of Ml 32 cells as observed 
during a scratch assay (Figure 9). 

2.8. IGF-1 Stimulation Induces Protein Expression of PRKCs in Ml 12 Osteosarcoma Cells 

It is well known that the insulin/insulin-like growth factor (IGF) pathway plays an important 
role in osteosarcoma tumor growth. Osteosarcoma tumors and cell lines express components of the 
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pathway (including insulin, IGF-I and -II, as well as pathway receptors) that are capable of autocrine 
signaling [44,45,51]. Additionally, inhibition of the pathway through various means is effective at 
inhibiting osteosarcoma growth in xenograft models [43,46-50]. IGF-I signaling is also known to lead 
to increased cellular levels of diacylglycerol, which can then activate PRKCs (and other proteins 
containing CI domains) [111-113], by a mechanism involving membrane tethering of PRKCe and 
conformational change [114]. Specifically PRKCs is activated following IGF-1 treatment in vascular 
smooth muscle cells, and may be involved in IGF-I-mediated proliferation and migration of these 
cells [112,113]. 



Figure 9. (a) PRKCe is not required for in vitro migration of highly metastatic Ml 32 
osteosarcoma cells, (b) Confirmation of PRKCe knockdown, (c) Quantification of migration 
by Ml 32 cells, (a) Upper panel: photographs of Ml 32 cells immediately after scratching, 
two days following transfection with PRKCs siRNA. Lower panel: matched photographs 
of Ml 32 cells two days following scratching and four days following transfection with 
PRKCe siRNA. (b) Western blot verification of PRKCe knock-down two, three and four 
days following transfection with siRNA. (c) Quantification of migrated distance by Ml 32 
cells transfected with PRKCe siRNA and appropriate controls. Values are average + standard 
error of three independent experiments. 
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The effect of IGF-1 stimulation on PRKCe protein in highly metastatic osteosarcoma cells was 
investigated. Highly metastatic Ml 12 osteosarcoma cells were serum starved for at least 24 h, followed by 
the addition of either fresh serum-free media or media containing 50 ng/mL IGF-1. It was observed 
that following IGF-1 stimulation, the protein expression of PRKCe increased in a time-dependent 
manner, with a peak occurring approximately 30 min following stimulation (Figure 10). 

Figure 10. Protein expression of PRKCs is induced by IGF-1 treatment. Western blots 
were performed on whole cell lysates of Ml 12 cells incubated in 1% FBS (low serum) 
medium or in 1% FBS medium supplemented with 50 ng/mL IGF-1 for the time periods 
described. Membranes were probed with anti-PRKCe antibody and anti-P-actin antibody as 
a loading control. The experiment was performed three independent times. 
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3. Discussion 

This work demonstrates that expression profiles of metastatic-at-diagnosis osteosarcomas (MD-OS) 
are quite distinct from expression profiles of localized-at-diagnosis osteosarcomas (LD-OS), possibly 
indicating a distinct disease etiology for the more severe MD-OS. Supervised network analysis 
discovered hundreds of networks that exhibited both differential activity and differential organization 
in the MD-OS expression profiles, at permissive cut-off levels. This indicates that either heterogeneous 
differences are observed in the MD-OS tumor group, or that changes to the transcriptome observed in 
MD-OS are potentially associated with many bystander events {i.e., aberrations in expression pattern 
not functionally relevant to osteosarcoma metastasis), or possibly both. 

Investigation of the cellular processes of networks within the permissive significant results showed 
differentially activated networks were strongly enriched for transport and translation networks, and 
slightly enriched for intracellular organization networks. Differentially organized networks were 
strongly enriched for metabolic networks, as well as intracellular organization and transport networks. 
Increased or over-active translation is known to support an aggressive phenotype of many cancers [115], 
and the analysis presented here implies it is a characteristic of MD-OS tumors as well. The enrichment 
of disorganized metabolic networks in MD-OS expression profiles may indicate that altered cellular 
metabolism plays a role in osteosarcoma metastasis, a notion supported by a recent study by Hua 
et al. [1 16]. By studying serum metabolite profiles of mice injected with OS cells, Hua et al observed a 
metabolic shift coincident with the onset of pulmonary metastasis [116]. 

Examination of the most stringently differentially activated and organized networks in MD-OS led 
to the selection of the PRKCs-RASGRP3-GNB2 network for follow-up characterization. This network 
was subsequently found to be differentially activated at the transcript level among a panel of in vitro 
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models of OS metastasis. PRKCe was also found to be more abundant in some highly metastatic OS 
cells at the protein level, and in some MD-OS tumors at the transcript level. Although PRKCs is 
known to support migration of other cell types [89,96,98,117,118], there was no evidence to support a 
role for PRKCe in the migration of Ml 32 osteosarcoma cells in this study. This may indicate that the 
aberrant expression of PRKCs and its network in vitro and in vivo is related to a bystander effect, or 
that the functional relevance of this aberrant expression has not yet been elucidated. The effect of high 
PRKCe expression may be related to the proliferation effects of the IGF-1 pathway, as this study 
provides evidence for the IGF-1 dependent induction of PRKCe protein expression. PRKCe may 
support other cellular processes entirely, as it is known to promote pro-metastatic phenotypes of many 
other cancers [82-87,91-103]. 

In addition to the PRKCs-RASGRP3-GNB2 and MYH9 network; this article describes several 
networks exhibiting significant differential activity, organization, or some combination of the two in 
MD-OS expression profiles. The integration of the results of this study with other datasets of osteosarcoma 
expression profiles, as they become available, will help to distinguish the drivers and genuine 
characteristics of metastatic osteosarcoma from the passengers, as would a limited high-throughput 
functional screen of the significant networks described in this study. 

4. Experimental Section 

4.1. Patient Follow -Up 

Overall survival data was available for all 46 patients in the group presenting without metastasis at 
diagnosis. Of these, 17 died of disease (DOD) with a median follow-up of 41 months (minimum 
follow-up = 6 month, maximum follow-up = 157, SD = 39.2 months), 28 are alive with no evidence of 
disease (ANED), with a median follow-up of 101 months (minimum follow-up = 35 month, maximum 
follow-up = 269, SD = 62.1 months) and one subject is alive with evidence of disease (AWED) with a 
follow-up of 12 months. 

Out of 17 patients presenting with metastases, overall survival data was available for only 14 
patients. Out these 14, 13 died of disease (DOD) with a median follow-up of 10 months (minimum 
follow-up = 1 month, maximum follow-up = 49, SD = 12.5 months) and the other patient is alive with 
no evidence of disease (ANED) with a follow-up of 178 months. 

4.2. Tumor Samples 

Primary high-grade intramedullary osteosarcoma tumors were selected for expression profiling by 
sarcoma pathologists on the basis of tumor homogeneity. The OS tumor cohort consisted of 63 tumors 
which were grouped into those that were localized (n = 46) or metastatic (n = 17) at the time of initial 
diagnosis. Samples were collected by open biopsies prior to administration of any chemotherapy, and 
stored in liquid nitrogen until time of RNA isolation. Total RNA was extracted using Trizol reagent 
(Invitrogen, Carlsbad, CA, USA). The amount and quality of RNA was assessed using both Ultrospec 
2100 pro (GE Healthcare Bio-Sciences, Piscataway, NJ, USA) and 1% agarose gels. 
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4.3. Gene Expression Profiling 

5 |ig of tumor and "reference" (pooled from cell lines) cDNA was indirectly labeled using 
aminoallyl nucleotide analogues with Cy3 and Cy5 fluorescent tags, respectively. The labeled cDNA 
was competitively hybridized to University Health Network 19 k cDNA arrays (UHN19k) containing 
18,981 "spots" (mapped to 8,998 known unique genes). This process was repeated with reciprocal 
fluorescent tagging. Data normalization, imputation (K10 Nearest Neighbours algorithm), and analysis 
were performed in collaboration with Drs. Shelley Bull, Dushanthi Pinnaduwage, and Robert Parkes. 
Supervised statistical analysis (random variance T-test) was performed by Robert Parkes using 
BRB-Array Tools software [119]. 

4.4. Unsupervised Hierarchical Clustering 

The most differentially expressed single genes within the expression profiling experiment were 
identified according to methods previously described [13,14]. In this study, a subset was examined that 
exhibited at least six-fold change in expression in at-least four tumors with a maximum p-vahxe 
of 0.001 (student's t-test). Unsupervised hierarchical clustering of the tumors according to their 
expression of these genes was performed with Partek Genomics Suite. 

4.5. Supervised Network Analysis 

The entire database of interactions for human genes and proteins was downloaded from the Pathway 
Commons website (as an adjacency list and the interactions were converted to Entrez GenelDs [15]. 
This dataset comprised physical and genetic interactions, as well as pathway and disease associations 
(which were either translated by Pathway Commons to binary interactions or were lost), from both 
curated and non-curated sources [15]. A subset of 5,855 genes was common to both this interaction 
database and the expression profiling experiment. There were 176,121 interactions among these genes, 
with six interactions being the median number per gene. 

4.5.1. Differentially Activated Networks 

Genetic networks demonstrating significant differential "activity" in MD-OS tumors were discovered 
in a manner similar, with some changes, to that previously described by Chuang et ah [11]. Briefly 
each network was restricted to those genes that met some cut-off of significant expression between 
localized and metastatic tumors (Table 5). 



Table 5. Cut-Off conditions for gene inclusion to "differential activity" analysis. 
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For each gene in the remaining network a "class difference score" (CDS) was calculated as the 
difference in median expression between classes (i.e., between metastatic tumors and localized 
tumors), and normalized to the variation within the localized samples [Equation (1)]. A "network 
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activity score" (NAS) was calculated by determining the average of the absolute CDS for each gene in 
the network [Equation (2)]. In these equations, gi, g2 ... gnare all genes in network J, which has Nq 
members. Ga and Gb are the median expression of gene G among localized (A) and metastatic (B) 
tumors. Sqa is the standard deviation of expression values of gene G among the localized tumors: 



cds 6 = 




(1) 



NASj = _ 

G=1l (2) 



To determine statistical significance, the NAS was compared to the corresponding NAS generated 
from 1,000 permutations of both sample and gene labels to determine two empirical ^-values 
[Equation (3)]. In these equations Jl . . . 0 - 1, 0, 0 + 1 ... J are all networks in the study, and i . . . t - 1, 
t, t + 1 ... Ni are the repetitions of gene and label permutations (in this study Ni = 1,000). 



r 0 i 
Pj = I f (NASj, NAS-Rjj) 



0 if NASj > NAS-Rji 
if NASj < NAS-Rji 
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The significance of the randomly generated NAS scores was also determined in a similar fashion 
[Equation (4)]; this was done for false-discovery rate (FDR) calculation. 
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Two FDRs were calculated for each network by determining the average number of 
randomly-generated networks with a p value equal to or lower than a nominal p value threshold [p^ in 
Equation (5) below], equal to the /?-value of the network being considered. The number of randomly 
generated ^-values falling below this nominal threshold was then divided by the number of real 
networks also falling below this nominal threshold [Equation (5)]. This was done for both /?- values 
(from sample and gene permutations) to yield two empirical FDR values for each network. The entire 
process was repeated for different cut-off conditions (Table 5), and was stopped when more stringent 
cutoffs failed to produce any additional significant networks: 
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4.5.2. Differentially Organized Network — Dynemo 

Genetic networks exhibiting significant differential organization in MD-OS samples relative to 
localized samples were discovered as previously described [12]. Briefly, for each network the Pearson 
Correlation Coefficient (PCC) was used to determine the overall correlation in gene expression 
between the network hub and each of its interactors in both the localized and metastatic samples 
[Equation (6)]. The difference in this hub-interactor correlation between localized and metastatic 
samples was then calculated for each interaction in the network [Equation (7)], and the average 
difference in correlation across the entire network was also calculated: "AvgAPCC" [Equation (8)]. In 
these equations Jl . . . 0 - 1, 0, 0 + 1 ... J are all the networks in the study. H is the hub (central node) 
of network J, which has Ng members, and gl, g2 . . . gn are all interactors of H, and therefore all other 
genes in network J. Sg and Sh are the standard deviations of gene G and hub H among the indicated 
tumor class, tl, t2 ... T are all the tumors in each class, and Nt is the total number of tumors in the 
localized (N T = 46) and metastatic (N T =17) classes: 

PCC = £ (GfG)(H t "H) 

1 1=1 (6) 

flPCCg = PCCg.M eTa5I3Tic - PCCg. LKB |j 2ec) (7) 
AvgAPCC = X I APCC sl 

E=gl - gn 



N, 



(8) 



The AvgAPCC value was compared to the corresponding scores generated from the same network 
following 1,000 permutations of the class labels to generate a non-parametric ^-value to assess the 
significance of the change in internal correlation of each network [Equation (9)]. Let i ... t - 1, t, 
t + 1 ... Ni be the number of repetitions of gene and label permutations (in this study Ni = 1,000). 

f 0 if AvgAPCCj > AvgAPCC-Rjj 
q = L 1 if AvgAPCCj < AvgAPCC-Rjj 

Pj = Z q (AvgAPCCj, AvgAPCC-Rjj) (9) 

j = J 

i = 1... N, 

N, 



4.5.3. Visualization of Network Results 



Two statistical confidence levels were investigated in this study, a permissive cut-off level to 
discover broad (or "global") trends among metastatic tumors, and a stringent cut-off to delineate 
high-confidence networks for follow-up. The permissive cut-off for differentially activated networks 
was chosen to be p's < 0.05 and FDR's < 0.2, and for the differentially organized levels the permissive 
cut-off was p < 0.001. These cut-offs were used to assess the significance of cellular process 
enrichment at the global level (Table 1) and for discovery of networks containing genes previously 
implicated in OS metastasis (Tables 3 and 4). Visualization of these global trends (Figure 1) was 
limited to computer processing power, and thus were further refined for the differentially activated 
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networks to be p's = 0.05, FDR's = 0.1 and NAS > 1.65, and for differentially organized networks to 
be p < 0.001 and |APCC| > 0.4. Visualization of the global trends was performed with the Enrichment 
Map plugin for Cytoscape [120]. The stringent cut-off levels for differentially activated networks were 
set to p's < 0.01, FDR's < 0.1, NAS > 1.65, and for differentially organized networks were set to p = 0, 
|APCC| > 0.4, and having at least seven interactors. Visualization of these high-confidence networks 
was conducted using Cytoscape [121] (Figures 3 and 4). 

4.5.4. Cellular Process Annotation 

Gene Ontology annotations (which relate genes to cellular processes) from the "Generic Slim" 
database were downloaded [122]. The database was further simplified to focus on interesting processes 
according to Table 6. Networks were assigned to processes by determining the most commonly 
occurring term among genes within the network. In this manner all 5,855 networks in the study could 
be assigned to a process (Figure 11a — the "study process composition", or Table 1, column 4). Since 
analysis of "differential activated" networks requires identification of significant subsets, all networks 
which were eventually found to be significantly differentially activated were re-annotated using only 
the genes within the network's significant subset (Table 1, column 2). This resulted in discordant 
annotations between the network and the significant subset for only 6% of all networks, and thus the 
annotations of the network subsets (Figure lib — "subsets process composition", or Table 1, column 2) 
and overall "study process composition" (Figure 11a, or Table 1, column 4) are overwhelmingly similar. 
As the methods of Taylor et al. do not identify significant subsets within networks, this consideration was 
not necessary for differentially organized networks. 

Table 6. Simplification of the gene ontology slim generic database. 

Original Terms Further Simplified Terms 

cell death, death death 



multicellular organismal development, embryonic 
development, anatomical structure morphogenesis 



development 



cell differentiation, differentiation differentiation 

regulation of gene expression, epigenetic epigenetics 

cell growth, growth growth 



cellular component organization, organelle organization, mitochondrion 
organization, cytoplasm organization 



intracellular organization 



metabolic process, cellular amino acid and derivative metabolic process, 

secondary metabolic process, lipid metabolic process, biosynthetic process, 

catabolic process, carbohydrate metabolic process, protein metabolic process, metabolism 

nucleobase nucleoside nucleotide and nucleic acid metabolic process, DNA 

metabolic process, generation of precursor metabolites and energy 

signal transduction, response to biotic stimulus, response to external stimulus, 

response to abiotic stimulus, cell-cell signaling, cell communication, response signaling 

to endogenous stimulus, cell recognition 

protein transport, transport transport 

regulation of biological process, biological process, behavior NaN 
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Figure 11. Distribution of cellular processes among networks and sub-networks investigated in 
this study, (a) "Study" Process Composition (b) "Subsets" Process Composition, (a) The 
Generic Slim database of gene-process associations was accessed from Gene Ontology and 
each network was annotated according to the most commonly occurring process term 
among genes within the network, (b) For all 497 networks found to be differentially 
activated in this study, the networks were re-annotated according to the most commonly 
occurring term among the significant sub-set (as the analysis of differential activity 
involves identification of a significant subset of genes within each network). 
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4.5.5. Cellular Process Enrichment 



Enriched cellular processes were determined separately for differentially activated and organized 
networks. Enriched processes were first identified by determining those processing comprising a 
greater proportion within the significant results than in the entire study (i.e., differentially activated: 
Table 1 column 7 was compared to Table 1 column 3, differentially organized: Table 1 column 10 was 
compared to Table 1 column 5). 

Significance of cellular process enrichment (/?-value) within the network results was determined 
using the hypergeometric distribution [Equation (10)], as described by Boyle et al. [123]. In this 
equation, N is the total number of networks in the study (5,855). M is the number of networks in the 
entire study annotated to a process of interest (differentially activated: Table 1 column 2, differentially 
organized: Table 1 column 4). i is the number of networks within the significant results annotated to 
the same process of interest (differentially activated: Table 1 column 6, differentially organized: Table 1, 
column 9), and n is the number of networks contained within that set of significant results (bottom of 
Table 1 — differentially activated: 497, differentially organized: 683). This operation was performed 
with Matlab using the "hygepdf command. 



P= 1 



/ M w N-M ) 
i n - i 



( N 

i 
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4.6. Cell Culture 



Hu09 and derivates (L06, L13, Ml 12 and M132) were grown in RPMI 1640 supplemented with 
10% fetal bovine serum (FBS) and 1% L-glutamine. L06 and L13 are cell lines derived from Hu09 
human OS cell line by limited dilution plating, and Ml 12 and Ml 32 were derived from the Hu09 line 
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by in vivo selection of pulmonary metastatic nodules [106,107]. L06 and LI 3 both form fewer lung 
metastases at a decreased incidence in mice following tail- vein injection compared to the intermediate 
Hu09 line, while Ml 12 and Ml 32 result in greater numbers of pulmonary nodules at a higher 
incidence [107]. The differences in metastatic propensity correlate to different survival lengths for 
mice injected with the various cell lines [107]. M8 cells were derived from MG63 human OS cells 
(hereafter referred to as MG63-P), and have decreased latency until pulmonary metastasis [109]. Both 
were grown in Dulbecco's Modified Eagle Medium supplemented with 10% FBS. LM7 cells were 
derived from SAOS2 human OS cells (hereafter referred to as SAOS2-P), and form pulmonary 
metastases at a greater incidence than SAOS2-P following tail-vein injection in mice [108]. LM7 and 
SAOS2-P were grown in McCoy's 5A supplemented with 5% FBS. RNA was harvested from all cell 
lines using Trizol reagent followed by phenol/chloroform extraction. 

4.7. Quantitative Reverse-Transcription Polymerase Chain Reaction (rt-PCR) 

RNA was collected from OS tumor samples as previously described [124], and cDNA was 
synthesized from both tumor RNA and cell line RNA according to the manufacturer's instructions 
(M-MLV Reverse Transcriptase — Invitrogen, Carlsbad, CA, USA). Quantitative rt-PCR was performed 
according to the manufacturer's instructions (Sybr Green — Applied Biosy stems, Life Technologies, 
Carlsbad, CA, USA) to quantify the abundance of target cDNA relative to that of a control gene, signal 
transducing adaptor molecule 2 (STAM2), using primer sequences according to Table 7. Dr. Dushanthi 
Pinnaduwage performed statistical analysis (Welch two-sample t-test on log2 transformed data) 
comparing expression levels of protein kinase C epsilon (PRKCe) transcript between localized and 
metastatic tumors. 

Table 7. Primer sequences used in this study. 

Gene Symbol Primer Pairs 

STAM2 Forward 5'-TGGATGACAGTGATGCCAATTG-3' 

Reverse 5'-CGCTGCCTCAGTCTCTATGT-3' 
PRKCe Forward 5 ' -C ACTGC AAGCTGGCTGACT-3 ' 

Reverse 5'-TGCAGGATCTCAGGAGCTATG-3' 
RASGRP3 Forward 5 ' -GG ATTTCTCTGGGGC AT A ATC- 3 ' 

Reverse 5'-AGGAGGTCTTTGCACTGTTTG-3' 
GNB2 Forward 5'-CTATCAAGCTGTGGGACGTG-3' 

Reverse 5'-GTAGCCGTTGGGGAAGAAAG-3' 
MYH9 Forward 5'-GCCTACAGGAGTATGATGCAAG-3' 

Reverse 5'-ACTGGATGACCTTCTTGGTGTT-3' 



4.8. Western Blots 

Cytosolic protein extracts were isolated from cell lines using NETN lysis buffer (150 mM NaCl, 
1 mM EDTA, 20 mM Tris pH 7.5, 0.5% NP40, 1 mM phenylmethylsulphonyl fluoride, and 1% each 
of protease inhibitor, phosphatase inhibitor I and phosphatase inhibitor II, all from Sigma-Aldrich, 
(St. Louis, MO, USA). Protein concentration was determined using the bicinchoninic acid (BCA) 
protein assay kit (Pierce, Thermo Scientific, Rockford, IL, USA). Proteins were separated using 10% 
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sodium dodecyl sulfate polyacrylamide gel electrophoresis (SDS-PAGE) and transferred to 
nitrocellulose membranes at 30 V, overnight at 4 °C. Membranes were blocked for 1 h with Tris 
buffered saline with 0.1% Tween-20 (TBS-T) supplemented with 5% fat- free milk. Primary antibody 
incubation was performed in TBS-T supplemented with 5% bovine serum albumin (BSA) or fat-free 
milk, according to the manufacturer's instructions. Secondary antibody incubation (donkey anti mouse, 
donkey anti rabbit) was performed at 1:5,000 concentration in TBS-T with 5% BSA for 50 min. 
Protein bands were visualized by chemiluminescence using ECL detection system (Amersham, GE 
Healthcare Bio-Sciences Corp, Piscataway, NJ, USA,). Primary antibodies used in this study are PKCe 
(1:1,000, Cell Signaling Technology, Danvers, MA, USA), and P-actin (1:5,000, Sigma-Aldrich, 
St. Louis, MO, USA). Analysis of PRKCe protein expression in the panel of OS cell lines was 
performed in triplicate. 

4.9. Knockdown ofPRKCs 

Two different siRNAs were purchased from Ambion (Invitrogen, Carlsbad, CA, USA) targeting 
PRKCe: "select si 1102" and "select si 1103" (hereafter referred to as PRKCs siRNA #2 and #3, 
respectively). Ambion "select negative control siRNA #1" was used as a negative control (hereafter 
referred to as "scramble"). Transfection reagents were purchased from Dharmacon Thermo Scientific 
(Rockford, IL, USA): "DharmaFECT 2" was used for the Ml 12 cell line, and "DharmaFECT 4" was 
used for the Ml 32 cell line. Transfections were performed in parallel with cell plating for experiments. 
20 |iM siRNA and the appropriate transfection reagent were each diluted (1.5:100 and 1:100, 
respectively) in Opti-MEM medium (Invitrogen, Carlsbad, CA, USA) and left to incubate for 5 min at 
room temperature (RT). The siRNA and transfection reagents were then mixed together and incubated 
for 10 min at RT. 200 |iL of the mixture was then applied to wells of 12-well plates, or 20 |iL was 
applied to wells of 96- well plates. 800 |iL (12-well plates) or 80 |iL (96- well plates) of cells at an 
appropriate concentration were then plated evenly in the wells. The cells were washed with phosphate- 
buffered saline (PBS) the following day and given fresh media. 

4.10. Scratch Assay 

Comparison of migration rate following knockdown of PRKCe in Ml 32 cells was performed by 
plating 1.2 x 10 5 cells in 12-well plates. All plates had grids drawn across them to allow repeated 
photographing of the same field. At least 12 fields were collected and analyzed for each sample on 
each day. Twenty four hours following plating, the cells were washed with PBS and given RPMI 1640 
supplemented with 1% each of FBS and L-glutamine ("low serum media"). On the second day 
following plating, the confluent cells were "scratched" using a 200 |iL pipette tip. The cells were washed 
twice with PBS and fresh low serum media was added before immediately imaging the cells using an 
inverted microscope and camera. The average distance that the confluent edge of cells had travelled 
into the wound was measured for each time point. This experiment was performed in triplicate. 
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4.77. IGF-1 Induction Assay 

To determine if IGF-1 is able to induce the expression of PRKCe, 2 x 10 5 Ml 12 or M132 cells were 
plated in 12- well plates. After two days the cells were washed twice with PBS and were given low 
serum RPMI media (1% FBS, 1% L-glutamine), with the exception of some cells which were retained 
in complete media as a control. The following day the remaining cells were given fresh low 
serum RPMI supplemented either with nothing ("no treatment" control) or with 50 ng/mL of IGF-1 
(Sigma- Aldrich, St. Louis, MO, USA). After the appropriate length of incubation with IGF-1, the cells 
(4 wells) were washed with cold PBS, scraped, and then lysed in order to harvest cytosolic protein as 
described above. Trizol reagent was added to a fifth well of each sample for the extraction of RNA and 
evaluation of mRNA abundance as described above. This experiment was repeated in triplicate. 

5. Conclusions 

Supervised network analysis was used to discover differentially activated and organized genetic 
networks in expression profiles of metastatic-at-diagnosis osteosarcomas (MD-OS) compared to 
localized-at-diagnosis osteosarcomas (LD-OS). The PRKCs-RASGRP3-GNB2 network was found to 
be differentially activated among MD-OS expression profiles and among in vitro models of OS 
metastasis. It was found that MD-OS tumors do not express significantly higher levels of PRKCe 
overall (t-test p = 0.1873), but they were more likely to exhibit high expression of PRKCs transcript, 
compared to the LD-OS tumors (five of fourteen MD-OS tumors had PRKCe expression greater than 
the maximum of the LD-OS tumors). 

This result is consistent with the expression pattern observed in the panel of OS cell lines, where 
PRKCs was found to be more abundant at the RNA level in some of the in vitro models of OS 
metastasis. Specifically, PRKCs was more abundant in the Ml 12 and M132 vs. Hu09, L06 and L13 
cell line model, but not in the LM7 vs. SAOS2 or M8 vs. MG63 models. The heterogeneity of PRKCs 
expression among MD-OS tumors may indicate heterogeneous networks are aberrant in 
MD-OS tumors, or that the PRKCs-RASGRP3-GNB2 network may be disrupted through alterations of 
other genes in the network. 

Despite reports by others of the involvement of PRKCs in cell migration [89,90,96,98,118], 
knockdown of PRKCs using siRNA was not found to affect migration of highly metastatic Ml 32 
osteosarcoma cells. The effect of PRKCs on invasion of osteosarcoma cells was not investigated in 
this study, and may be a fruitful avenue of further research, as PRKCs is known to be involved in 
invasion of other cell systems in vitro, as well as to be involved in various other pro-metastatic 
pathways [125,126]. The absence of pro-metastatic effect of PRKCs protein in migration assays 
indicates that either the aberrant expression of the PRKCs-RASGRP3-GNB2 network may be related 
to a bystander effect, or that the pro-metastatic phenotype of PRKCs over-expression remains to 
be elucidated. 

PRKCs protein expression in highly metastatic Ml 12 cells was found to be induced by IGF-1 
stimulation, and this may indicate that PRKCs is involved in IGF-1 dependent pathways. Beyond the 
PRKCs-RASGRP3-GNB2network, this article describes many aberrantly activated and organized 
networks among expression profiles of MD-OS tumors. A systematic functional screen of these 
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networks, and determination of the predictive accuracy of these network expression patterns in 
independent datasets, would help differentiate the bystanders from the drivers of OS metastasis. 
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